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Abstract. This is an example of C++ code of Clifford algebra calculations with the GiNaC computer algebra system. 
This code makes both symbolic and numeric computations. It was used to produce illustrations for paper ]l4[ 

Described features of GiNaC are already available at PyGiNaC jsj and due to course should propagate into other 
software like GNU Octave pf and gTybalt MM which use GiNaC library as their back-end. 
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1. Introduction 

This example of Clifford algebras calculations uses GiNaC library Q], which includes a support for generic Clifford 
algebra starting from version 1.3.0. Both symbolic and numeric calculation are possible and can be blended with 
other functions of GiNaC. Described features of GiNaC are already available at PyGiNaC || and due to course should 
propagate into other software like GNU Octave [|| and gTybalt |L8| which use GiNaC library as their back-end. 

We bind our C++-code with documentation using noweb fl6|| within the literate programming concept JDJ. Our 
program makes output of two types: some results are typed on screen for information only and the majority of 
calculated data are stored in files which are lately incorporate by MetaPost |ll[ to produce PostScript graphics for 
the paper |lj, |l2). Since this code can be treated as software we are pleased to acknowledge that it is subject to 
GNU General Public License ffij. 

GiNaC allows to use a generic Clifford algebra, i.e. 2 n dimensional algebra with generators efc satisfying the 
identities CiCj + eye, = B(i,j) + B(j,i) for some (metric) B(i,j), which may be non-symmetric ||, ^ and contain 
symbolic entries. Such generators are created by the function 

ex clifford_unit(const ex & mu, const ex & metr, unsigned char rl = 0, bool anticommuting — false); 
where mu should be varidx class object indexing the generators, an index mu with a numeric value may be of type 
idx as well. Parameter metr defines the metric B(i, j) and can be represented by a square matrix, tensormetric or 
indexed class object, optional parameter rl allows to distinguish different Clifford algebras (which will commute each 
other). The last optional parameter anticommuting defines if the anticommuting assumption (i.e. e^e^ + e^e^ = 0) 
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will be used for contraction of Clifford units. If the metric is supplied by a matrix object, then the value of 
anticommuting is calculated automatically and the supplied one will be ignored. One can overcome this by giving 
metric through matrix wrapped into an indexed object. 

Note that the call cliff ord_unit(mu, minkmetric()) creates something very close to dirac_gamma(mu) , although 
dirac_gamma have more efficient simplification mechanism. The method Clifford:: get_metric() returns metric defining 
this Clifford number. The method clifford::is_anticommuting() returns the anticommuting property of a unit. 

If the matrix B(i,j) is in fact symmetric you may prefer to create the Clifford algebra units with a call like that 

ex e = cliff or dju,nit{mu, indexed(5, sysymmQ, i, j)); 
since this may yield some further automatic simplifications. Again, for a metric defined through a matrix such a 
symmetry is detected automatically. 

Individual generators of a Clifford algebra can be accessed in several ways. For example 

{ 

varidx nu(symbol("nu"), 4); 

realsymbol s("s"); 

exM = diag_matrix(lst(l, -1, 0, s)); 



ex 


e = 


- cliff ord_unit(nu, 


ex 


eO 


= e.subs{nu = 


0) 


ex 


el 


= e.subs{nu = 


1) 


ex 


e2 


= e.subs{nu = 


2) 


ex 


e3 


= e.subs(nu = 


3) 



} 

will produce four generators of a Clifford algebra with properties e\ = 1, ef = — 1, e| = and e| = s. 
A similar effect can be achieved from the function 

ex lst_to-clifford(const ex & v, const ex & mu, const ex & metr, unsigned char rl = 0, 
bool anticommuting = false); 

ex lst-to-clifford(canst ex & v, const ex & e); 
which converts a list or vector v= (vq, v±, . . . , v n ) into the Clifford number vgeo + v±ei + ■ • • + v n e n with e& directly 
supplied in the second form of the procedure. In the first form the Clifford unit is generated by cliff or d_unit(mu, 
metr, rl, anticommuting). The previous code may be rewritten with help of lst_to-clifford{) as follows 

{ 

varidx nM(symbol("nu"), 4); 
realsymbol s("s"); 
exM= diag-matrix(\st(l, -1, 0, s)); 
ex eO = lst_to_clifford(lst(l, 0, 0, 0), nu, M) 
ex el = lst_to_clifford\\st{0, I, 0, 0), nu, M) 
ex e"2 = lst_to_clifford\\st{0, 0, 1, 0), nu, M) 
ex e,3 = lst_to_clifford{lst(0, 0, 0, 1), nu, M) 

y 

There is the inverse function 

1st cliff ord-to Js£(const ex & e, const ex & c, bool algebraic=true) ; 
which took an expression e and tries to find such a list v= (v ,vi, . . . ,v n ) that e = v cq + V\C\ + • • • + v n c n with 
respect to given Clifford units c and none of Vk contains the Clifford units c (of course, this may be impossible). 
This function can use an algebraic method (default) or a symbolic one. In algebraic method Vk are calculated as 
(ecfc + Cke)/pow(ck,2)- If pow(ck,2) is zero or is not numeric for some k then the method will be automatically 
changed to symbolic. The same effect is obtained by the assignment {algebraic— false) in the procedure call. 

There are several functions for (anti-)automorphisms of Clifford algebras: 

ex clifford_prime(const ex & e) 

inline ex clifford_star(const ex & e) { return e.conjugateQ; } 

inline ex clifford_bar(const ex & e) { return cliff ord_prime{e.conjugate()); } 
The automorphism of a Clifford algebra cliff ord-prime() simply changes signs of all Clifford units in the expression. 
The reversion of a Clifford algebra cliff vrdstarQ coincides with conjugate^) method and effectively reverses the 
order of Clifford units in any product. Finally the main anti-automorphism of a Clifford algebra clifford_bar{) is the 
composition of two previous, i.e. makes the reversion and changes signs of all Clifford units in a product. Names for 
this functions corresponds to notations e', e* and e used in Clifford algebra textbooks Q, ||. 
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The function 

ex clifford_norm(const ex & e) ; 

calculates the norm of Clifford number from the expression ||e|| 2 = ee. The inverse of a Clifford expression is returned 
by the function 

ex clifford-inverse(const ex & e) ; 

which calculates it as e" 1 = e/ ||e|| 2 . If |je|j = then an exception is raised. 

If a Clifford number happens to be a factor of dirac_ONE() then we can convert it to a "real" (non-Clifford) 
expression by the function 

ex remove-dirac-ONE(const ex & e) ; 
The function canonicalize_clifford() works for a generic Clifford algebra in a similar way as for Dirac gammas. 

The last provided function is 

ex clifford_moebius-map(const ex & a, const ex & b, const ex & c, const ex & d, const ex & v, const ex & G, 

unsigned char rl = 0, bool anticommuting = false); 

ex clifford_moebius-map(const ex & M, const ex & v, const ex & G, unsigned char rl = 0, 

bool anticommuting = false); 
It takes a list or vector v and makes the Mobius (conformal or linear-fractional) transformation Q 



v i ► (av + b) (cv + d) 1 defined by the matrix M 



a b 
c d 



The matrix may be given in two different forms — as one entity or by its four elements. The last parameter G define 
the metric of the surrounding (pseudo-)Euclidean space. This can be an indexed object, tensormetric, matrix or 
a Clifford unit, in the later case the optional parameters rl and anticommuting are ignored even if supplied. The 
returned value of this function is a list of components of the resulting vector. 
Finally the function 

char clifford_maxJabel(const ex & e, bool ignore-ONE = false); 
can detect a presence of Clifford objects in the expression e: if such objects are found it returns the maximal 
representationJabel of them, otherwise -1. The optional parameter ignore-ONE indicates if dirac_ONE objects 
should be ignored during the search. 

I^Tj5]X output for Clifford units looks like \clif f ord [1] {e}"{{\nu}}}, where 1 is the representationJabel and \nu 
is the index of the corresponding unit. This provides a flexible typesetting with a suitable defintion of the \clif f ord 
command. For example, the definition 

\newcommand{\clif f ord} [1] [] {} 
typesets all Clifford units identically, while the alternative definition 

\newcommand{\cliff ord} [2] []{\ifcase #1 #2\or \tilde{#2} \or \breve{#2} \fi} 
prints units with representation_label=0 as e, with representation_label=l as e and with representation Jabel=2 as e. 

2. Main procedure 

Here is the main procedure, which has a very straightforward structure. This and next initialisation section is 
pretty standard. The first usage of GiNaC for Clifford algebras is in Section^. 

<*!>= S> 

(Includes 0) 

(Definitions 

(Global items jrt|) 

int main (int argc, char**ar<7«) 

{ 

(C++ variables declaration |sa|) 
(CiNaC variables declaration] 
(Pictures tuning |a|) 



(Parabola parameters 19c } 



Defines: 

main, never used. 
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Now we run a cycle over the three possible type of metric in two dimensional space (i.e. elliptic, parabolic and 
hyperbolic). For each space we initialise the corresponding Clifford units, symbolically calculate various types of 
Mobius transforms as well as vector fields for three subgroups of SL 2 (R). 

«|£]> 

for (metric — elliptic; metric < hyperbolic; metric++) { 

cout <C endl <^endl <C "Metric is: " <C metric-name[metric] <C " . " <C enrfZ; 
(Initialise Clifford numbers ^ 



(Calculation of Moebius transformations 10c > 



(Calculation of vector fields lit) 



Uses elliptic |Ha| hyperbolic and metric |7b| . 

Then we run a cycle for three subgroups of SL 2 (K) (i.e. A, iV, if). For all possible combinations of those with 
metric from the surrounding cycle in the previous chunk we 

(1) build orbits of the subgroups and their transversal curves; 

(2) two types of the Cayley transform images of all above curves; 

(3) check some formulae in the paper; 

We draw all pictures by substitution of numeric values into the symbolic results obtained in the above chunks. 

(*!>+= «@E3> 

for (subgroup — subgroup_A; subgroup < subgroup_K; subgroup++) { 
/* iteration over su bgro ups A, N and K */ 



(Drawing arrows 12b) 
(Building orbits 12c) 
(Building transverses 13d ) 

} 

} 

Uses subgroup |7b| , subgroup_A |Ha| and subgroup_K 

Finally we draw eight frames which illustrates the continuous transformation of the future part of the light cone into 
the its past part |l^, Figure [|. 

<0 



(Build future-past transition 14c) 

} 



3. Auxiliary matter 
Some standard inclusions, but do not forget GiNaC library! 

(Includes 0)= (|) 
^include <ginac/ginac .h> // At least ver . 1.4.0! 
^include <cmath> 

using namespace std\ 
using namespace GiNaC] 
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3.1. Defines. Some constants are defined here for a better readability of the code. 

(Definitions |j| = (§) |i]> 

/ / Defined constants 
#define elliptic 
#define parabolic 1 
#define hyperbolic 2 
#define subgroupjl 
#define subgroup_N 1 
#define subgroup_K 2 
#define grey . 6 



18b 



Defines: 

elliptic, used in c hunk s [ia[ [k| , and 
grey, used in chunk |l2b| . 

hyperbolic, used in chunk s nA 6aj |1C|, 14—16, and|l8dl 



parabolic, used in chunks [LObl and 1 7—19 . 
subgroup_A, used in chunks 11 , O, |l5c] , 16c , and ^(j. 



subgroup_K, used in chunks 
subgroup_N, used in chunk H6 



and 11-19. 



Some macro definitions which we use to make more compact code. They initialise variables, open and close cm 
description in the MetaPost file. Here is initialisation of a new curve. 

(Definitions Q+= (§) |c]o 

#define init_coord(X) upos [X] =0; \ 



vpos[X] = 
udir[X] = 1 
vdir[X] = 



fprintf(fileout[X], "draw ") 



Defines: 

init_coord, used in chunks ^ and 15s. 
Uses f ileout and upos [ja[ 

This part is used to close a curve output in cases the end is reached or curves passes the infinity. 

(Definitions |a|}+= (|) |d|o 

#define close.curve (X) fprintf (f ileout [X] , " (a°/„+5 . 3f ,b°/„+5 . 3f ) *u withcolor 7.5 . 3f *7.s ; \n" , \ 

upos[X], vpos[X], color_grade, color_name[subgroup\) 
#define put_draw(X) fprintf (f ileout [X] , "\ndraw ") 



Defines: 

close_curve, used in chunks ^d, |l4b , and |l5a 



put_draw, use d in chunks I 



Uses color_name gaj, fileout] 



sub 



ubgroup [7b| , u |8aJ and upos [saj 

Here is the common part of code which is used for outputs a segment to files. 

(Definitions 0>+= (§) <i|c] |a|> 

#define put_point(X) if (inversion) \ 

fprintf(ftleout[X\, "(a7.+5.3f ,b7.+5.3f )*u. . .", upos[X] , vpos[X]); \ 
else if (direct) \ 

fprintf(fileout[X], "(a7.+5.3f ,b7.+5.3f )*u{(7.+5.3f ,7.+5.3f )}. .", upos[X] , vpos[X], udir[X], vdir[X\); \ 
else \ 

fprintf(fileout[X], " (a7.+5 . 3f ,b7.+5 . 3f ) *u{ (7.+5 . 3f ,7.+5 . 3f ) } . . .", upos[X] , vpos[X], transjuf, trans.vf); 
^define renew_curve (Y) close_curve(Y) ; \ 

put_draw( Y) 



Defines: 

put_point, used in chunk 
renew_curve, used in chunk 

Uses close_curve ^c| f ileout 



put_draw ^c| u ^a|, and upos 
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We should make a rough check that the curve is still in the bounded area, if it cross infinity then such line should 
be discontinued and started from a new. Our bound are few times bigger that the real picture, the excellent cutting 
within the desired limits is done by MetaPost itself with the clip currentpicture to . . . ; command. 

Besides some outer margins we put different types of bound depending from the nature of objects: sometimes it 
is limited to the upper half plane, sometimes to hyperbolic unit disk. The necessity of such checks in the hyperbolic 



case is explained in |l2j, § 2.5 1 



(Definitions Q+= (|) <[k] |l]> 

#define if _in_limits (X) if ( (abs(ujres .to_double() ) <= ulim) \ 

A (abs(v_res.to_double()) < vlim) \ 
A ((metric ^ hyperbolic) \ 

V inversion \ 

V ( -icayley A (v_res.is-positive() V vjres.is-zero())) \ 

V ( cayley A ^ex_to<num.eric>(-pow(u-res,2)+pow(v-res,2)-1.001).is-positive()))) { \ 

iy?os[X] = u-res.to-double(); \ 
vpos[X] = v_res.to_double(); \ 
ex Vect = dV[subgroup][X\.subs(lst(x = u_res, y= v_res)); \ 
udir[X] — ex Jo<numeric> (Vect. op(0)).to-double(); \ 
vdir[X] = ex_to<numeric>(Vect.op(l)).to_double(); \ 
putjpoint(X) ; \ 
} else { \ 

renew Lcurve(X) ; \ 



Defines: 

if _in_limits, used in chunks 15-17. 
Uses hyperbolic ^a|, metric (fb| numeric |sd| , put_point |jd| , renew_curve subgroup 0, ulim Q and upos pa| . 

Then a curve is going through infinity we catch the exception, close the corresponding draw statement of MetaPost 
and start a new one from the next point. 

(Definitions Q+= (§) [kjo 

#define catchJiandle (X) cerr << "*** Got problem: " << p.whatO << endl ; \ 

renew _curve(X) 

Defines: 

catch_handle, used in chunks 15—17. 
Uses renew_curve ^d| 

Extracting of numerical values out of Moebius transformations 

(Definitions Q+= (|) f^t> 

^define get_components u_res = ex_to<numeric> (res . op(0) . evalf () ) ; \ 

V-res — ex-to<numeric>(res.op(l).evalf()) 



Defines: 

get_components, used in chunks 15—17. 
Uses numeric pdL 
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To make an accurate drawing we calculate the direction of a transverse line out of symbolic vector fields calculation 
done before. 

(Definitions Q+= (§) <^c| [l7a|t> 

#define transverse_dir (X) if (! direct) { \ 

trans-uf = ex- to<numeric> ( trans-dir_sub[X] .op(0).subs( ajnode) . evalf()) . to_ doubleQ ; \ 

trans-vf = ex_to<numeric>(trans-dir_sub[X\.op(l).subs(a_node).evalf()).to_double(); \ 

if (trans_uf = INFINITY) { trans _uf = 1; trans _vf = 0; > \ 

else if {trans_uf= -INFINITY) { trans.uf = -1; trans_vf = 0; > \ 

else if {trans_vf= INFINITY) { trans_uf =0; trans.vf = 1; } \ 

else if (transjvf= -INFINITY) { trans_uf=0: trans.vf = -1; } \ 

else if (abs(trans-uf)+ abs(trans-vf) > 100) { \ 

double r= sqrt(trans-uj * trans _uf + transjuj '* trans_vf); \ 

trans_uf -i-= r; trans_vf -j-= r; \ 
>\ 

> 



Defines: 

transverse di r. used in chunks hq and Il7l 
Uses numeric pd| . 

Some global variables which is convenient to use for the openfileQ function below. 

(Global items 0) = (|) 0> 

int subgroup, metric; / / Subgroup iterator and sign of the metric of the space 
numeric signum; 

char sgroup[] = "ANK", metric_name[) = " eph" ; // Names used for a readable output 



6-10. and 13-19. 



lib 



19c , and 20 



Defines: 

metric, used in chunks la, 
sgroup, used in chunks 7c 
subgroup, used in chunks 4—9 and 11-20. 
Uses numeric |jd| . 

Procedure openfileQ is used for opening numerous data files with predefined template of name. 

(Global items Q+= (|) «([b] 

FILE * openfile(const char *F) { 

char filename[} = " cayley-t-k-e.d", fempZ[] = "cayley-t- , / . Is- /, ls.d"; 
char *5 filename— filename, * Stempl=templ; 
strcat(strcpy(Stempl, F), "-'/,. Is-'/, . Is .d"); 

sprintf(Sfilename, Stempl, lksgroup[subgroup], &zmetric_name[metric])\ 
return fopen(Sfilename, "w"); 

} 



Defines 

openf ile, 
Uses metr 



1 e. used in chunks [12] and 13d . 

ic |7bl sgroup [7bJ, and subgroup |7b|. 



8 VLADIMIR V. KISIL 



1st February 2008 



3.2. Variables. First we define variables from the standard C++ classes. 

(C++ variables declaration |a|} = (^) 
FILE *fileout[S\; /* files to pass results to \MetaPost\ */ 

static char * color _name\\~{ "hyp", "par", "ell", "white"}, // A, N, K subgroup colours 
=f=/ormM/a[]=-["\nDistance to center is:", "\nDirectrice is:", 

"\nDif f erence to foci is:"}; 
bool direct = true, // Is it orbit or transverse? 

cayley — false, // Is it the Cayley transform image? 

inversion = false; //Is it future-past inversion? 
double u, v, upos[3], vpos[3], udir[3], vdir[3], vval = 0, // coordinates of point, vector, etc 

color-grade, focaLJ[2] = {0, 0}, 

trans juf=l, trans-vf—0; 



Defines: 



12b 



15a 



color_name, used in chunks jx and 
flleout, used in ch unks j^, jll |lj , and 
u, used in chunks ^, |Ljb| , an* 
upos, used in chunksWand I 
v, used in chunks pi [12, and 
Uses subgroup [7b] . 

Then other variables of GiNaC types are defined as well. They are needed for numeric and symbolic calculations 

(CiNaC variables declaration |t]) = (|) |c]> 

varidx mt(symbol("nu", "\\nu"), 2), mM(symbol("mu", "\\mu"), 2), 

psi(symbol("psi", "\\psi"),2), xi(symbol("xi", "\\xi"), 2); 
realsymbol x("x"), y("y"), t("t"), // for symbolic calculations 

a("a"), 6("b"), c("c"), // parameters of the parabola v = au 2 + bu + c 

tr_u("U"), tr_v("V"); / / Vector ot the tranverse direction 
1st ajnode, soln[2], a-trans; 

Uses v [ja| , 

(CiNaC variables declaration ^k])+= (|) <^ ^t> 

matrix M(2, 2), // The metric of the vector space 

C(2, 2), Cl(2, 2), CI(2, 2), C1I(2, 2), // Two versions of the Cayley transform 
T(2, 2), 17(2, 2), // The map from first to second Cayley transform 

Jacob[3] [3]={matrix(2, 2)}, // Jacobian of the Moebius transformation 
trans-dir[3] [3]={matrix(2,l)}, / / Components of transverse direction 

irans_<izr_sM6[3]={matrix(2,l)}; // Components of transverse direction substituted by a subgroup 

Uses metric and subgroup [fb| 

(CiNaC variables declaration gt])+= (§|) «|sc] |<]> 

numeric u_res, v_res, //coordinates of the Moebius transform 

up[3][2] = {0, 0, 0, 0, 0, 0}, // saved values of coordinates of the parabola 

vp[3][2] = {0, 0, 0, 0, 0, 0}; 
const numeric halffl, 2); 

Defines: 

numeric, used in chunks ^, 0, Q |l2b| , and 18—20. 

(CiNaC variables declaration 0)+= (|) <0 

ex res, e, eO, el, 

Moebius[3][5], dV[3}[3}, 
II indexed by [subgroup], [type](=direct, Cayley-operator, Cayley 1-op, C-point, Cl-p) 
rfc?y[3], Curv[3], I) indexed by [type](=direct, Cayley-operator, Cayleyl-op) 
focal, p, 

focaLl, focaLu, focaLv; / / parameters of parabola given by v = au 2 + bu + c 



Uses subgroup |7b| and v 
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Here is the set of constants which allows to fine tune MetaPost output depending from the type of subgroup and 
metric used. The quality of pictures will significantly depend from the number of points chosen for iterations: too 
much lines will mess up the picture, very few makes it incomplete or insensitive to the singular regions. These 
numbers should be different for different combinations of subgroup and metric. 

(Pictures tuning ^ = (|) 
const int m/irmis[3] [3] = {10, 20, 30, // indexed by [subgroup], [metric] 
10,10,19, 
10, 10, 10}, 

/sfeps[3][3] = {15, 15, 20, // indexed by [subgroup], [metric] 
15, 10, 20, 
12, 15, 15}; 
float ulim = 25, vlim = 25, 

flimits[3] [3] = {2.0, 2.0, 4.0, // indexed by [subgroup], [metric] 
10.0, 4.0, 4.0, 
0.5, 0.5, 0.5}, 

vpoints[3][10] = {0, 1.0^-8, 1.0-=-4, l.CK-2, 1.0, 2.0, 3.0, 5.0, 8.0, 16.0,// [metric] [point] 
0, l.O-r-8, 1.0H-4, I.Oh-2, 1, 2.0, 3.0, 6.0, 10.0, 20.0, 
0,1.0^8, 1.0-=-4, 1.0H-2, 1.0, 2.0, 3.0, 5.0, 10.0, 100}; 



Defines: 



Ik 



ulim, used in chunks |6aj and 
vilimits, used in chunks 13—16. 
Uses metric and subgroup J?b| 

Initialise the set of coordinates for a cycle 

(Initialisation of coordinates = (13a 14a) 

init_coord(0); 
init-Coord{\); 
put_draw(2); 

Uses init_coord ^ and put_draw |Hc]. 

4. Symbolic Clifford Algebra Calculations 

This section finally starts to deal with Clifford algebras. We try to make all possible calculations symbolically 
delaying the numeric substitution to the latest stage. This produces a faster code as well. 

4.1. Initialisation of Clifford Numbers. We initialise Clifford numbers first. Alternative ways to define e, eO 
and el are indicated in comments. 

(Initialise Clifford numbers (0@) ^> 

signum — numeric(mefric-l); / / the value of e\ 
M=-l, 0, 
0, signum; 

//e = clifford_unit(mu, indexed(M, symmetric2(),xi, psi)); 

e = cliff ord_unit(mu, M); 

eO = e.subs{mu = 0); 

el = e.subs{mu = 1); 

// e0= lst_to_clifford(lst(1.0, 0), mu, M); 

II el = lst_to_clifford\\st{Q, 1.0), mu, M); 

Uses metric [7b| and numeric jjcl . 

Now we define matrices used for definition of the alternative Cayley transforms JlJ, [l^] . 

(Initialise Clifford numbers 0)+= (£a 14c) <[k] [l0a| t> 

T = dirac-ONEQ, eO, / / Transformation to the alternative Cayley map 

f 1 ei 

eO, dirac-ONE(); //is given by I 1 

y Cr 1 

TI= dirac-ONEQ, -eO, // The inverse of T 

-eO, dirac_ONE(); // is given by f e ^ 
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A form of matrix for the Cayley transform depends from the type of metric. 

(Initialise Clifford numbers 0)+= (^e 14c) «|id| 

switch (metric) { 
case elliptic: 
case hyperbolic: 

1 -e 2 



10b > 



C = dirac_ONEQ, -el, / / First Cayley transform . ^ 
signurmel, dirac_ONEQ; 

1 e2 



01 = dirac.ONEQ, el, II The inverse of O . 

V — a e2 1 

-signurmel, dirac_ONEQ; 

CI = C.mul(T); / / Second Cayley transform 

C1I= TI.mul{CI); // The inverse of CI 

break; 

Uses elliptic |iaj hyperbolic and metric |7h| . 

In the parabolic case there are two different (elliptic and hyperbolic) types of the Cayley transform |lj, [l^] . 

(Initialise Clifford numbers |jk]}+= ([la 14c) < 10a 

case parabolic: 

C = dirac_ONEQ, -el*half, / / First (elliptic) Cayley transform for the parabolic case 

-el*half, dirac.ONEQ; 
CI= dirac.ONEQ, ehhalf // The inverse of C 

eUhalf dirac.ONEQ; 

CI = dirac.ONEQ, -el* half, / / Second (hyperbolic) Cayley transform for the parabolic case 

eUhalf dirac-ONEQ; 
C1I= dirac-ONEQ, eHhalf // The inverse of CI 

-el*half dirac.ONEQ; 
break; 



Uses elliptic |ia| hyperbolic and parabolic ^a| 



4.2. Mobius Transformations. We calculate all Moebius transformations along the orbits as well as two their 
Cayley transforms only once in a symbolic way. Their usage will be made through the GiNaC substitution mechanism. 
First, we define matrices Exp.A, Exp.N, Exp.K Jl5|, VI. 1] related to the Iwasawa decomposition of SL 2 (M.) 
§111.1]. 



(Calculation of Moebius transformations 10c ) = 
matrix Exp.A(2, 2), Exp.N(2, 2), Exp.K{2, 2); 

Exp.A 











exp(t) * dirac-ONEQ, 0, // Matrix 

0, exp(-t) * dirac.ONEQ; 

Exp.N = dirac.ONEQ, t * eO, / / Matrix ( 

0, dirac-ONEQ; 

Exp.K = cos(t) * dirac.ONEQ, sin(t) * eO, / / Matrix 

sin(t) * eO, cos(t) * dirac.ONEQ; 
ex Exp[3] = {Exp.A, Exp.N, Exp.K}; 
matrix E(2, 2); 







= exp 



exp 

1 

C 



( cos t — sin t 
\ sin t cos t 



exp 



-1 

1 
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Here we symbolically calculate the related Mobius transformations, as well as its two images under Cayley transform 
C and CI for operators and finally two Cayley images for points. 



10c 



(Calculation of Moebius transformations !0c)+= (£; 
for {subgroup = subgroup_A; subgroup < subgroup_K; subgroup++) { 
try{ 

Moebius[subgroup][0] = cliff ord_moebius_map(Exp[subgroup], 1st (a;, y), M); 
/* Cayley transforms of operators */ 

Moebius[subgroup] [1] = clifford_moebius_map(canonicalize_clifford( 

C.mul(ex_to<matrix>(Exp[subgroup]).mul(CTj)), lst(x, y), M); 

Moebius[subgroup] [2] = cliffordjmoebiusjmap(canonicalize_clifford( 

Cl.mul(ex_to<matrix>(Exp[subgroup]).mul(ClTj)), lst(x, y), M); 

I '* Cayley transforms of points */ 
Moebius[subgroup] [3] = clifford_moebius_map(C.mul(ex_to<matrix>(Exp[subgroup])), lst(a;, y), M); 

Moebius[subgroup] [4] = clifford_moebius_map(canonicalize_clifford( 
Cl.mul(ex_to<matrix>(Exp[subgroup]))), 1st (a;, y), M); 
} catch (exception &cp) { 

cerr <C "*** Got problem in vector fields: " <C p.whatQ <C endl; 



Uses catch 



12a 



16b 



17b 



17c, subgroup |7b|, subgroup_A Egl and subgroup_K [ja] 



4.3. Symbolic Calculations of the Vector Fields. We c alcu late sym bolic expressions for the vector fields of 
three subgroups A, N and K, which are stated in |l^, Lemmas 2A and 2^2 1 and shown on Jl2|, Figures |and§. The 



formula for a derived representation p(X) of a vector field X used here is § VI. 1] 

d 
Jt 



p(X) = ^-p(e tx ) | t=0 . 



Results of calculations are directed to stdout and will be used for a better drawing of orbits, see QA . 

We calculate Jacobi&n of the Mobius transformations in order to suply to MetaPost the tangents of the transverse 
lines. 

(Calculation of vector fields Ub)= (0) [l2a| i> 



try { 

cout < "Vect field \t Direct \t\t In Cayley \t\t In Cayleyl" < endl; 
for (subgroup — subgroup-A; subgroup < subgroup-K; subgroup++) { 
for (int cayley = 0; cayley < 3 ; cayley++) { 

1st Moeb = ex_to<\st>(Moebius[subgroup][cayley! cayley+2:0\); 
dV[subgroup] [cayley] — Moebius[subgroup][cayley].diff(t).subs(t = 0); 
Jacob[subgroup] [cayley] = matrix(2,2, lst(Moeb.op(Q).diff(x), Moeb.op(0).diff(y), 

Moeb.op(l)-dtff(x), Moeb.op^.diffiy))); 
1 1 Transformation of a direction by a Jacobian 
for (int i=0; i < 2; i++) 

trans _dir[subgroup][cayley] — Jacob[subgroup][cayley].mul(matri^c(2, 1, lst(£r_«, tr_v))); 

} 

cout <C " d" <C sgroup[subgroup] <C " is:\t" <C dV[subgroup][0] 
<C ";\t" <C dV[subgroup][l] <C " ; \t " < dV[subgroup] [2] < end^; 

> 

Uses sgroup [7^ subgroup [?bj subgroup_A ^a|, and subgroup_K ^a|. 
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We also calculate curvature of the orbits using the formula [§, § 5.1(20)] 

\xy - yx\ 



k = 



(x 2 + y 2 y' 2 



(Calculation of vector fields llb)+= (|4aj) 
for (int cay ley = 0; cay ley < 3 ; cayley++) { 

ddV[cayley] = Moebius[subgroup_K\[cayley].diff(t, 2).subs(t = 0); 
ex du = ex-to<lst>(dV[subgroup_K][cayley]).op(0); 
ex dv — ex-to<lst>(dV[subgroup-K\[cayley]).op(l); 
ex ddu = ex_to<lst>(ddV[cayley\).op(0): 
ex ddv = ex_to<lst>(ddV[cayley\).op(l); 

Curv[cayley] — normal((ddu * dv - du * ddv)^rpow(du*du+dv*dv, 1.5)); 

} 

cout <C "Curvature of K-orbits on the v-axis: " 

<C nor7na/(Cury[0].su6s(x = 0)) <C eniiZ; 

} catch (exception Szp) { 

cerr -C "*** Got problem in vector fields: " <C p.what() <C endl; 

y 



lib 



Defines: 

catch, used in chunks Lla[. 15b, 19a, and po) 
Uses subgroup_K [5a| and v 

5. Numeric Calculations with Clifford Algebras 
Numeric calculations are done in to fashions: 



(1) Through a substitution of numeric values to symbols in some previous symbolic results, subsection 5.1 and 



(2) Direct calculations with numeric GiNaC classes, subsection 5.3. 
The first example of substituion approach is the drawing of vector fields. Three vector fields are drawn by arrows 
into a MetaPost file. 



(Drawing arrows 12b >= ( jibj ) 



fileout[0] — openfile(" arrows"); // open MetaPost file 
color_grade = grey; 

for (int k = -10; k < 10; k++) j / cycle over horizontal dir 
for (int j = 0; j < 11; { // cycle over vertical dir 
u = fc-^3.0; / / Calculate coordinates of the point 
v = £-3.0; 

u_res = ex-to<numeric>(dV[subgroup][0].subs(lst(x = u, y = v)).op(Q)); j j Get numeric 
vjres = ex_to<mimeric>(dV[subgroup\[0].subs(lst(x = u, y = v)).op(l)); 

fprintf(fileout[0] , "myarrow ( (a%+5 . 3f , b°/ +5 . 3f ) , (a'/.+5 . 3f %+5 . 3f *s , b°/.+5 . 3f °/.+5 . 3f *s) ) " , 

u, v, u, v,-res.to-double(), v, V-res.to-doubleQ); 
fprintf(fileout[0], " withcolor °/ 5 . 3f *°/ s;\n", color_grade 1 color _name[Z\)\ 

} 

fclose(fileout[0]); 

Uses color_name ^a|, f ileout grey numeric |id|, openf ile subgroup (7b| u and v |iaj 

5.1. Numeric Calculations of Orbits and Transverses. For any of three possibility e\ = —1, 0, 1 and three 
possible subgroups (A, N, K) orbits are constructed. First we open output MetaPost files. 
(Building orbits |l2c} = @ [l3a| t> 



direct = true; 
fileout[0] = openfile(" orbit"); 

fileout[l] — openfile(" cayley"); // Cayley transform of the orbits 

fileout[2] — openfile(" cayl-&"); / / Alternative Cayley transform of the orbits 

Uses f ileout ^ and openf ile [rc|. 
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This chunk runs iterations over the different orbits, which are initiated by the point vi. 



(Building orbits |l2c}+= @ 
for (int vi = 0; vi < vilimits[subgroup] [metric] ; vi++) { // iterator over orbits 
color-grade = l.2*vi-^-vilimits[subgroup][metric]; 
if (subgroup = subgroup_K) 

cout <C formula[metric] ; 
(Initialisation of coordinates ^ 
(Nodes iteration s |l3b| ) 
(Close all curves |i4b|) 
(Check parabolas 19c) 

} 



12c 



(Closing all files 13c) 



Uses metric subgroup subgroup_K ^5a|, and villmits ^)a| 

Each orbit is processed by iteration over the "time" parameter j on the orbit. For each node on an orbit an entry 
is put into appropriate MetaPost file, formulae from papers are numerically checked and all Cayley transforms are 
produced. 

(Nodes iterations |l3b| ) = 

for (int j = -fsteps[subgroup] [metric] ; j < fsteps[subgroup][metric]; j++ ) { 



(13a) 



float / = flimits[subgroup 



(Generating one entry 15c) 



[metric] *j~fsteps[subgroup] [metric] ; / / the angle of rotation 



(Check formulas in the paper 18a) 



(Producing Cayley transform of the orbit 17b) 



} 

Uses metric [7^ and subgroup (?b| 

Closing all files when finishing drawing orbits or transverses 



(Closing all files |l3c| ) 
fclose(fileout[0] ) ; 
fclose(fileout[l]); 
fclose(fileout[2]); 



(13a 14a) 



Uses f ileout pa| . 

5.2. Building of Transverses. Construction of transverses to the orbits follows the same structure as for orbits 
themselves with the changed order of iterations over time parameter and orbit origins. We again start from opening 
of the corresponding files. 

(Building transverses 13c) = (0) 14a > 

direct = false; 

fileout[0] — openfile(" orbit-t"); 

fileout[l] — openfile(" cayley-t"); // Cayley transform of transverses 
fileout[2] = open/ife("cayl-a-t"); // Alternative Cayley transform of transverses 
color --grade = 1.2; 



(Define transverse directions |l6c| ) ; 



Uses f ileout pA and openf ile 
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Thus chunk performs iterations over the transverse lines. 

14a (Building transverses |l3d| )+= 

for (int j = -fsteps[subgroup][metric]; j < fsteps[subgroup][metric\; j++ ) { 

float /= flimits[subgroup\[metric]*j^-fsteps[subgroup\[metric];/ / the angle of rotation 
(Initialisation of coordinates ^ 

for (int vi — 0; vi < vilimits[subgroup][metric]; vi++) { // iterator over orbits 
vval = vpoints[metric][vi]; 
(Generating one entry 15c) 



(Producing Cayley transform of the orbit 17b) 



(Close all curves |l4b| ) 

} 

(Closing all files 13c) 



Uses metric [/b| subgroup [fb| and vilimits 

All MetaPost draw statements should be closed at the end of run. 



(13a 14a) 



141 (Close all curves |i4 
close-curve(0) 
close_curve{l) 
close-Curve(2) 

Uses close_curve 

5.3. Future-to-Past Transformations. We finish our code by producing a set of frames of transformations from 
future to past of light cone. First we create a necessary hyperbolic setup and use subgroup_K to fill in parts of the 
light cone. 

14c] (Build future-past transition |l4c|) = (^) 



15a > 



(Build future-past transition 14c ) = 

const int curves =15, nodes = 40, frames = 8; 
const float expscale — 1.3, nodescale = 4.0, 

rad[] = {1.0-=-5, 1.0-r4, l-=-3.5, 1.0-=-3, l-=-2.5, 1.0-=-2, l-=-1.5, 1, 1.5, 2, 2.5, 3, 3.5, 4, 4.5, 5} 
char name\\ — "future-past-OO.d"; 
char* S = name; 
ulim = 8.5; 
vlim = 8.5; 



metric — hyperbolic; 
subgroup — subgroup_K; 
direct = true; 

inversion = true; //to cut around unit circle. 
(Initialise Clifford numbers 0) 

ex Fut = clijford_moebius-map(dirac-ONE() 1 -a * el, 

a * el, diracONEQ, lst(x, y), M); 



Defines: 

curves, used in chunk |15 
exp_scale, used in chunk 
name, never used 
used in chunk |15a 



Uses hyperbolic ^a, metric (?b| subgroup [7b], subgroup_K [5a|, and ulim |9a} 
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Then we proceed with iteration over the frames. 

(Build future-past transition 14c )+= (0) <i |l4c| 

for (int j = 0; j < frames; j++) { // the power of transformation 
sprintf(S, "future-past-°/o.2d.d", j); 
fileout[0] = fopen(S, "w"); 

for (int k — 0; k < curves; k++) { // the number of curve 
color_grade = k^-frames; 
init_coord(0); 

(Iteration over a curve L5b> 
close-curve(0); 

} 

fclose(fileout[0]); 

y 

cout <C endl; 



llc- 



Uses close_curve curves 14c, fileout gaj init_coord gbj and S 

A curve is created by rotation in hyperbolic metric and then a Mobius transformations corresponding to the frame 
number is applied. 

(Iteration over a curve 



I5q) = 

for (int I — -nodes~2; I < nodes+2; 1++) j j the node number 
try { // There is a chance of singularity! 

float angl = ((j > 0) ? exp(douhle(j^ exp_scale-3)) : 0); 
res = Fut.subs(lst( a = angl, x = rad[k]* cosh(l-^-node_scale) , y 
get_components; 
if_injimits(0); 
} catch (exception 8zp) { 
catch_handle(Q); 



(15a) 



rad[k]* sinh(l-^ nodes cale) ) ) ; 



Uses catch 12a 16b 17b 17c, catchjiandle |jb[ exp_scale 14c, get_components |jc], and if _in_limits [fe]. 



5.4. Single Node Calculation. This is a common portion of code for building orbits and transverses. A single 
entry into MetaPost file is calculated and written. Calculations depend for three possible values of the subgroup. For 
each of them we create a special list ajnode of substitutions for the already symbolically calculated Moebius[] [] . 
For subgroup A the points are distributed evenly on the unit circle. 

(Generating one entry 15c >= (13h| 14a) l5do 

switch (subgroup) { 
case subgroup-A: 

vval = l.Q*vi^r(vilimits[subgroup\[metric]-l); 

if (metric = hyperbolic) //we need a double set of value for negatives as well 
vval *= 2; 

ajnode = lst(£ = /, x = cos(Pi*vval), y = sin(Pi*vval)); 
break; 



Uses hyperbolic metric [7bJ subgroup [7^, subgroup.A [5a|, and vilimits |9a| 

For the subgroups K and N points are distributed evenly on the vertical axis. 



(Generating one entry 15c >+= 
case subgroup_K: 

vval — vpoints[metric][vi]; 
ajnode = 1st ( t = (/ * Pi) , x ■ 
break; 



(13b 14a) < 



15c 



16a > 



0, y = vval); 



Uses metric and subgroup_K 
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For the subgroup N we need additional points. 



(13b 14a) < 



15d 



16b > 



(Generating one entry 15c >+= 
case subgroup-N: 

if (metric = hyperbolic) {//we need a double set of value for negatives as well 
vval = ( ((vi - vilimits[subgroup] [metric] -j-2) < ) ? -1 : 1) 
*vpoints[metric][abs(vi - vilimits[subgroup\\metric] -j-2)]; 

} else 

vval = vpoints[metric][vi]; 
a_node = lst(t = /, x = 0, y = vval); 
break; 

> 



Uses hyperbolic ^5a|, metric [?b| subgroup [7^, subgroup_N |3a|, and vilimits |ia| 

And now using values stored above in ajnode we do the actual calculation through substitution and write the node 
into the MetaPost file. 



(Generating one entry 15c )+= 
try{ 

res — Moebius[subgroup\[0].subs(a-node); 
transverse_dir(0) ; 
get-components; 
if_in_limits(0) ; 
} catch (exception &p) { 
catch_handle(0); 

y 



(13b 14a) < 



16a 



Defines: 

catch, used in chunks [lla , 15b, 19a| , and |2p[ 
Uses catch_handle ^ get.components g if .in.limits g subgroup @ and transverse.dir 

We define the tangents to the transverse lines here. 



(Define transverse directions 16c >e 
switch (subgroup) { 
case subgroup-A: 

a_trans — 1st (tr_u = -y, tr_v 

break; 
case subgroup_N: 
case subgroup-K: 

a_trans — 1st (tr_u = 0, tr_v : 

break; 

> 



(13d) 



16d> 



x); 



i); 



Uses subgroup [7b] , subgroup_A |5a| subgroup_K and subgroup_N |ia 

And here again comes evaluation through substitution. 



(Define transverse directions 16c )+= 
for (int cal=Q; cal < 3; cal++) 
trans-dir_sub[cal\ = matrix(2, 1, 

1st ( trans_dir[subgroup] [ ca\ .op(0).subs( a_trans) . normalQ, 
trans_dir[subgroup][cal\.op(l).subs(a_trans).normal())); 

Uses subgroup |7b| . 



(13d) < 



16c 
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5.5. Cayley Transforms of Images. We will need a calculation parameters of parabola into the Cayley transform 
images. This checks the statement |l2| Lemma 2.17]. The calculation is done by linear equation solver IsolveQ from 
GiNaC. 

(Definitions Q+= (|) <pa| [l9b|> 

^define calc_par Jocal (X) if (direct && (metric == parabolic) \ 

A (subgroup ^ subgroup_K)) { \ 
wp[2][X] = ujres; \ 
vp[2][X] = vjres; \ 

if (J = 1) i \ 
1st eqns, vars ; \ 
vars = a, 6, c; \ 
eqns = a*pow(up[0][X], 2) + b*up[0][X] + c = up[0][X], \ 

a*pow(up[l][X}, 2) + b*up{l][X] + c = vp[l][X], \ 

a*pow(up[2][X}, 2) + b*up[2][X] + c = vp[2}[X]; \ 
soln[X] — ex_to<lst>(lsolve(eqns, vars)); \ 

> \ 
/* After calculation is made we store previous values for the next round. */\ 
up[0][X] = up[l][X]; \ 
vp[0][X] = vp[l][X}; \ 
up[l][X] = up[2][X\; \ 
vp[l][X] = vp[2][X}; \ 



Defines: 

calc_par_f ocal, used in chunk [nj. 
Uses metric [7b| parabolic subgroup frb[ and subgroup_K [sa| . 

We produce two versions of the Cayley transforms for each node of the orbit or transverses lines. This done by 



(13b 14a) 



17c > 



simple substitution of ajnode into Moebius[] [3,4] symbolically calculated in subsection 4.2 

(Producing Cayley transform of the orbit |l7b[ > = 
cayley = true; 

try { // There is a chance of singularity! 

res ~ Moebius[subgroup][3].subs(a_node); 
transverse-dir(l); 

get_components; 

if_inJimits(l); 

calc-par_focal(Q) ; 
} catch (exception hp) { 

catch_handle(\); 

> 

Defines: 

catch, used in chu pks H 1 a 
Uses calc_par_f ocal 



15b 



19e 



20 



17a, catch_handle pb[ get_components [ic| if _in_limits subgroup [7b|, and transverse_dir 



For second type of the Cayley transforms we perform an extra run similar to the above. 



(Producing Cayley transform of the orbit 17b )+= 
try { 

res = Moebius[subgroup][A].subs(a_node); 
transverse-dir(2)\ 
get_components\ 
if_inJimits(2); 
calc-par_focal( 1 ) ; 
} catch (exception hp) { 
catch_handle(2); 

y 

cayley = false; 



(131 14a) 



17b 



Defines: 



catch, used in chunks Lla 



15b 



19^ , and |o[ 



Uses calc_par_f ocal 17a, catchjiandle ^b| get_components |jc| if _in_limits subgroup [7b|, and transverse_dir 
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5.6. Numeric Check of Formulae. Here is a numeric check of few formulas in the paper about radius and focal 
length of sections. We calculate focal properties for three types of orbits (circles, parabolas and hyperbolas) of the 
subgroup _K. 



(Check formulas in the paper 18a >e 



(13b) 



18b > 



if ((j ^ -fsteps[subgroup][metric]) A (j ^ fsteps[subgroup][metric]) // End points are weird! 
A (subgroup = subgroup_K) A (vval ^ 0)) {// only for that values 
try { 

switch (metric) { // depends from the type of metric 

Uses metric [7b| subgroup [?b| and subgroup_K |ia| 



The values are calculated for each node on the orbit as follows, see |12| , Lemma 2.2 \. For elliptic orbits of sungroup_K: 
a circle with the centre at (0, (v + v~ 1 )/2) and the radius (v — v~ 1 )/2. 



(13b) < 



18a 



18c > 



(Check formulas in the paper 18a )+= 
case elliptic: 

focaLj[l] = ex_to<mxmev'ic>(pow(u-res*ujres 

+ pow(v-res-(vval+l-^-vval)^-2, 2) , 0.5)) .to-double(); 

break; 

Uses elliptic |5a| and numeric gg . 

For parabolic orbits of sungroup_K: a parabola with the focus at (0, (v + v~ l )/2) and focal length v^ 1 /2. 



18d> 



(Check formulas in the paper 18a )+= 
case parabolic: 

focaLj[l] = ex_to<mimeric>(pow(u-res*U-res 

+ pow(v-res-(vval+l^-vval-^4), 2) , 0.5)-V-res).to-double(): 

break; 

Uses numeric |8d| and parabolic |ia] 

For hyperbolic orbits of sungroup_K: a hyperbola with the upper focus located at (0, /) with: 

for < v < 1: and 



/ = 



P 



P- 



1. 



Or — I, for V > 1. 



and has the focal distance between focuses 2p. 



(13b) < 



18c 



19a > 



(Check formulas in the paper 18a >+= 
case hyperbolic: 

p = (vval*vval+l)^-vval-^-pow(2,0.5); 

focal = ((vval<l) ? p -pow(p*p J ^2-l,0.5) : p+pow(p*p~2-l,0.5)); 
focaLj[l] — ex_to<numeric>(pow(u-res*U-res + pow(v_res-focal, 2) , 0.5) 
-pow(u-res*U-res + pow(v_res-focal+2*p, 2) , 0.h)).to-double(): 

break; 

} 

Uses hyperbolic and numeric ^d|. 
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If the obtained value is reasonably close to the previous one then = sign is printed to stdout, otherwise the new value 
is printed. This produce lines similar to the following: 

Distance to center is: 3.938====================== 

Remark 5.1. Note that all check are passed smoothly (see Appendix [A|), however in the hyperbolic case there is "V" 
shape of switch from positive values to negative and back (with the same absolute value) like this: 
Difference to foci is: 2.000====== -2.000============== 2.000====== 

This demonstrates the non-invariance of the upper half plane in the hyperbolic case as explained in [fL2l 



2.5]. 



(13b) < 



18d 



(Check formulas in the paper 18a >+= 

if ((abs(focaLf[l] - focaLf[0}) < 0.001) 

V (abs(focaLj[l] - focaLJ[0]) < OMl*(abs(focaLftL])+abs(focaLj{0})))) 
cout < " = "; 
else 

printf( " 7,5 . 3f " , focaLfil] ) ; 
focaLM = focalj[l]; 
} catch (exception Szp) { 

cerr <C "*** Got problem in formulas: " <C p.whatQ <C endl; 

y 



Uses catch 



12a 



16k 



17b 



17c 



The last check we make is about some properties of Cayley transform in parabolic case. All parameters of parabolic 



orbits were calculated in Subsection 5.5, now we check properties listed in jl2|, Lemma 2.17] 
(Definitions |ia|)+= 



^define output Jocal (X) ex_to<numeric>(focal_u.subs(soln[X] ) .evalf ()) .to_double() , \ 

ex-to< numer ic> (focaLv. subs(soln[X] ) . evalf( )).tO-double(), \ 
ea;_£o<numeric> (focal J. subs(soln[X\ ) . evalfQ). to_double() 



Defines: 

output_f ocal. used in chunk 19d 
Uses numeric [sd|. 

Here is expressions for focal length focaLl and vertex (focaLu, focaLv) of a parabola given by its equation v 
au 2 + bu + c. 



(Parabola parameters 19c >= (g) 
focaLl — l-=-(4*a); // focal length 
focaLu = b+(2*a); // u comp of focus 
focaLv = c-pow(b-^(2*a), 2); // v comp of focus 

Uses u |ja| and v paj 

Properties of parabolas are printed to stdout in the form: 

Parab (A/ 7/ 0.368); vert=( 1.140, -2.299); 1= 0.2500; second vert=(-l . 140, -2.299); l=-0.2500 



The two vertexes correspond to two Cayley transformations P e defined by C and P h defined by CI, see [|12| § 2.6] 
(Check parabolas 19d)= (Ei3) ^> 



if ((metric = parabolic) A (subgroup ^ subgroup_K)) 
try { 

printfC AnParab ('/, . ls/°/,2d/°/, 5.3f); vert=('/. 5.3f, % 5.3f); l=°/„ 5.4f", 

Szsgroup[subgroup], vi, vval, outpuLfocal(0)); 
printf("; second vert=(7„ 5.3f, °L 5.3f); l=°/„ 5.4f", outputjocal(l)); 
Uses metric output_f ocal 19t, parabolic sgroup subgroup [7b|, and subgroup_K [se]. 
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In the case of subgroup_A an additional line 
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Check vertices: -1 and -1 

is printed. It co nfirm s that vertexes of the orbits under the Cayley transform do belong to the parabolas v — ±v 2 — 1, 
as stated in [|l| 



(13a) < 



19d 



(Check parabolas l9d}+= 

if (subgroup = subgroup_A) 

cout <C "\nCheck vertices: " 

•C ex_to<mimeric>(focaLv.subs(soln[0\) + pow(foeal_u.subs(soln[0]), 2).evalj()).to_double() 
< " and " 

<C ex_to<nuxneric>(focal_v.subs(soln[l]) + pow(focal_u.subs(soln[l\), 2) .evalfQ) .to_double(); 
} catch (exception Szp) { 

printf( "\nParab ('/. . ls/ / 2d/°/ 5 . 3f ) is a straight line", 

&zsgroup[subgroup], vi, vval); 

y 



Uses catch 



12a 



16k 



17b 



17c, numeric |jd[ sgroup |7^, subgroup |7b| and subgroup_A ^a|. 
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6. How to Get the Code 

(1) Get the tyr^X source of this paper fl3|| from the |arXiv.org . 

(2) Run the source through P/TjtX. Three new files (noweb, C++ and MetaPost sources) will be created in the 
current directory. 

(3) Use it on your own risk under the GNU General Public License ]l0[ . 

Appendix A. Textual Output of the Program 
Here is the complete textual output of the program. 



Metric is: e. 

Vect field Direct In Cayley In Cayleyl 

dA is: 2*x,2*y; -2*y*x, l-y~2+x~2; -l-y~2+x~2,2*y*x 

dN is: 1,0; l/2-y+l/2*y~2-l/2*x~2,-y*x+x; -y*x-y , l/2+x-l/2*y~2+l/2*x~2 
dK is: l-y~2+x~2,2*y*x; -2*y,2*x; -2*y,2*x 
Curvature of K-orbits on the v-axis: (l+y~4-2*y~2) " (-1 . 5) * (-2*y~5-2*y+4*y~3) 

Distance to center is: 

Distance to center is: 3.938====================== 

Distance to center is: 1.875====================== 

Distance to center is: 0.750====================== 

Distance to center is: 0.000====================== 

Distance to center is: 0.750====================== 

Distance to center is: 1.333====================== 

Distance to center is: 2.400====================== 

Distance to center is: 3.938====================== 

Distance to center is: 7.969====================== 

Metric is: p. 

Vect field Direct In Cayley In Cayleyl 

dA is: 2*x,2*y; 2*x,2+2*y+2*x~2; 2*x, 2+2*y-2*x~2 
dN is: 1,0; l,2*x; l,-2*x 

dK is: l+x~2,2*y*x; l+x~2,2*y*x+4*x; l+x~2,2*y*x 
Curvature of K-orbits on the v-axis: -2*y 

Parab (A/ 0/ 0.000); vert=( 0.000, -1.000); 1= 0.2500; second vert=( 0.000, -1.000); l=-0.2500 
Check vertices: -1 and -1 

Parab (A/ 1/ 0.053); vert=( 0.083, -1.007); 1= 0.2500; second vert= (-0.083, -1.007); l=-0.2500 
Check vertices: -1 and -1 

Parab (A/ 2/ 0.105); vert=( 0.172, -1.029); 1= 0.2500; second vert=(-0. 172, -1.029); l=-0.2500 
Check vertices: -1 and -1 

Parab (A/ 3/ 0.158); vert=( 0.271, -1.073); 1= 0.2500; second vert=(-0. 271, -1.073); l=-0.2500 
Check vertices: -1 and -1 

Parab (A/ 4/ 0.211); vert=( 0.389, -1.151); 1= 0.2500; second vert=(-0.389, -1.151); l=-0.2500 
Check vertices: -1 and -1 

Parab (A/ 5/ 0.263); vert=( 0.543, -1.295); 1= 0.2500; second vert= (-0.543, -1.295); l=-0.2500 
Check vertices: -1 and -1 

Parab (A/ 6/ 0.316); vert=( 0.765, -1.586); 1= 0.2500; second vert=(-0.765, -1.586); l=-0.2500 
Check vertices: -1 and -1 

Parab (A/ 7/ 0.368); vert=( 1.140, -2.299); 1= 0.2500; second vert=(-l . 140, -2.299); l=-0.2500 
Check vertices: -1 and -1 

Parab (A/ 8/ 0.421); vert=( 1.974, -4.898); 1= 0.2500; second vert=(-l . 974, -4.898); l=-0.2500 
Check vertices: -1 and -1 

Parab (A/ 9/ 0.474); vert=( 6.034, -37.410); 1= 0.2500; second vert=(-6 . 034, -37.410); l=-0.2500 
Check vertices: -1 and -1 

Parab (A/10/ 0.526); vert=(-6.034, -37.410); 1= 0.2500; second vert=( 6.034, -37.410); l=-0.2500 
Check vertices: -1 and -1 

Parab (A/11/ 0.579); vert=(-l . 974, -4.898); 1= 0.2500; second vert=( 1.974, -4.898); l=-0.2500 
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Check vertic 
Parab (A/ 12/ 
Check vertic 
Parab (A/ 13/ 
Check vertic 
Parab (A/ 14/ 
Check vertic 
Parab (A/ 15/ 
Check vertic 
Parab (A/ 16/ 
Check vertic 
Parab (A/17/ 
Check vertic 
Parab (A/ 18/ 
Check vertic 
Parab (A/ 19/ 
Check vertic 
Parab (N/ 0/ 
Parab (N/ 1/ 
Parab (N/ 2/ 
Parab (N/ 3/ 
Parab (N/ 4/ 
Parab (N/ 5/ 
Parab (N/ 6/ 
Parab (N/ 7/ 
Parab (N/ 8/ 
Parab (N/ 9/ 
Directrice i 
Directrice 
Directrice 
Directrice 
Directrice 
Directrice 
Directrice 
Directrice 
Directrice 
Directrice 



es: -1 and -1 

0.632); vert=(-1.140, 
es: -1 and -1 

0.684); vert= (-0.765, 
es : -1 and -1 

0.737); vert= (-0.543, 
es : -1 and -1 

0.789); vert= (-0.389, 
es : -1 and -1 

0.842); vert=(-0.271, 
es: -1 and -1 

0.895); vert=(-0. 172, 
es: -1 and -1 

. 947) ; vert= (-0 . 083 , 
es: -1 and -1 

1.000) ; vert=(-0.000, 
es: -1 and -1 
0.000) 
0.125) 
0.250) 
0.500) 
1.000) 
2.000) 
3.000) 
6.000) 
10.000) 
20.000) 



vert=( 0.000, - 
vert=( 0.000, - 
vert=( 0.000, - 
vert=( 0.000, - 
vert=( 0.000, 
vert=( 0.000, 
vert=( 0.000, 
vert=( 0.000, 
; vert=( 0.000, 
; vert=( 0.000, 



-2.299) ; 1= 0.2500; 

-1.586) ; 1= 0.2500; 

-1.295) ; 1= 0.2500; 

-1.151) ; 1= 0.2500; 

-1.073) ; 1= 0.2500; 

-1.029) ; 1= 0.2500; 

-1.007) ; 1= 0.2500; 

-1.000) ; 1= 0.2500; 



•1.000) 
0.875) 
•0.750) 
0.500) 
0.000) 
1.000) 
2.000) 
5.000) 

9.000) ; 

19.000) 



second vert=( 
second vert=( 
second vert=( 
second vert=( 
second vert=( 
second vert=( 
second vert=( 
second vert=( 



1.140, 
0.765, 
0.543, 
0.389, 
0.271, 



-2.299) 
-1.586) 
-1.295) 
-1.151); 1=- 
-1.073); 1=- 



0.172, -1.029); 1= 
0.083, -1.007); 1= 
0.000, -1.000); 1= 



0.2500 
0.2500 
0.2500 
0.2500 
0.2500 
0.2500 
0.2500 
0.2500 



1= 





2500 


second 


vert= 


( 





000, 


-1 


1= 





2500 


second 


vert= 


( 





000, 


-0 


1= 





2500 


second 


vert= 


( 





000, 


-0 


1= 





2500 


second 


vert= 


( 





000, 


-0 


1= 





2500 


second 


vert= 


( 





000, 





1= 





2500 


second 


vert= 


( 





000, 


1 


1= 





2500 


second 


vert= 


( 





000, 


2 


1= 





2500 


second 


vert= 


( 





000, 


5 



1= 0.2500; second vert=( 0.000, 
1= 0.2500; second vert=( 0.00C 



.000) 


1= 


-0. 


2500 


.875) 


1= 


-0. 


2500 


.750) 


1= 


-0. 


2500 


.500) 


1= 


-0. 


2500 


.000) 


1= 


-0. 


2500 


.000) 


1= 


-0. 


2500 


.000) 


1= 


-0. 


2500 


.000) 


1= 


-0. 


2500 


9.000); 1 


=-c 


.2500 


19.000) ; 


1= 


-0.2500 



s : 

s: 1.875=== 

s: 0.750=== 

s: 0.000=== 

s: -0.750== 

s: -1.875== 

s: -2.917== 

s: -5.958== 

s: -9.975== 

s: -19.987= 



Metric is: h. 

Vect field Direct In Cayley In Cayleyl 
dA is: 2*x,2*y; -2*y*x, l-y~2-x~2; 2*y,2*x 

dN is: 1,0; l/2-y+l/2*y~2+l/2*x~2,y*x-x; l/2-l/16*(16*y-16*x)*x+l/2*y~2-l/2*x-2, l/2-l/16*(16*y-16*x)*j 
dK is: l+y~2+x~2,2*y*x; l+y~2+x~2,2*y*x; l+y~2+x~2,2*y*x 
Curvature of K-orbits on the v-axis: (-2*y~5-2*y-4*y~3) * (l+y~4+2*y~2) " (-1 . 5) 



Difference 


to 


foci 


is 


Difference 


to 


foci 


is 


Difference 


to 


foci 


is 


Difference 


to 


foci 


is 


Difference 


to 


foci 


is 


Difference 


to 


foci 


is 


Difference 


to 


foci 


is 


Difference 


to 


foci 


is 


Difference 


to 


foci 


is 


Difference 


to 


foci 


is 



8 . 125 -8 . 125========================== 8 . 125 

4.250= -4.250======================== 4.250= 

2.500=== -2.500==================== 2.500=== 

2.000====== -2.000============== 2.000====== 

2.500========= -2.500======== 2.500========= 

3.333========== -3.333====== 3.333========== 

5.200============ -5.200== 5.200============ 

10.100============= -10.100 10.100=========== 

100.010============= -100.010 100.010======== 



Appendix B. A Sample of Graphics Generated by the Program 

A sample of graphics produced by the program and post-processed with MetaPostis shown on Figure [l| Some 
more examples can be found in . 
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Figure 1. The elliptic, parabolic and hyperbolic unit disks. 



Appendix C. Index of Identifiers 



calc_pa r_f ocal: |17a| , |17b| , [17c 

" |l^1i6U Ir7b 



11a 



catch: 

catchJiandle 
close_curve: |5d. I5d 



e: 16 bl 
: 5dJ 



14c 



15bUl6bUl7b] 



color_na me: jcj |8a 



141 



12b 



15a 



17c 



15a 



elliptic: [la lja, [l0a| ^0b| , [l8b 



exp_scale: |14cj |15 
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12bl 12c. 13c, 13d, 15a 




16b, 17b. 17c 



fileout: |5b|, [ik], |5d 
get_components: |6d |15 
grey: Q ph 
hyperbolic: |la|, jj, |6a| [l0a| , |l0b| p^ , |l5c| ^6aj |l8d 
if_in_limits: |pj]l5b|, |l6bj |l7b|, [l7cf 
init_coord: [5b|, |9b|, [l5a| 
main: 



numeric: 



metric : ga| , |f3a|, |7b|, [7c], g fk], [TOa], |l3^, [l3t], ||a], |l4c], |l5^ [D^, [l6aj [l7a|, ||aj |9d 
14cj 

,R R pbl , pb] pdj |gg , |2 



6a| 

openf ile: |7_c 
output jfocal 

parabolic: |a|, |10b|, |17gj , |18c| , |19 
put_draw: 
put_point 
renew.curve 
S: pic] , 
sgroup 
s ubgroup: 
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